Introduction to Spatial Analysis
1 Rasters
1.1 Downloading raster climate data
Worldclim database http://worldclim.org/version2
# using the raster package
bio10 <- getData("worldclim",var="bio",res=10, path="./data/")
## Group of raster layers
class(bio10)## [1] "RasterStack"
## attr(,"package")
## [1] "raster"## class : RasterLayer
## dimensions : 900, 2160, 1944000 (nrow, ncol, ncell)
## resolution : 0.1666667, 0.1666667 (x, y)
## extent : -180, 180, -60, 90 (xmin, xmax, ymin, ymax)
## coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
## data source : in memory
## names : bio1
## values : -26.9, 31.4 (min, max)
1.2 What is a raster?
“Gridded” data that represents areas on the Earth’s surface as pixels.
1.3 Exploring a raster
1.3.1 Spatial resolution
The size of each cell in meters
1.3.2 Spatial extent
X, Y coordinates of the corners of the raster in the geographic space
1.3.3 Coordinate Reference System & Projection Information
The projection of a dataset refers to the way the data are “flattened” in into a 2D space. The coordinate system references the x and y coordinate space that is associated with the projection.
“If you have the same dataset saved in two different projections, these two files won’t line up correctly”
1.4 Plotting raster
1.4.1 Changing breaks and colours
## $breaks
## [1] -40 -20 0 20 40
##
## $counts
## [1] 11693 181932 203695 187201
##
## $density
## [1] 0.001000221 0.015562486 0.017424096 0.016013197
##
## $mids
## [1] -30 -10 10 30
##
## $xname
## [1] "v"
##
## $equidist
## [1] TRUE
##
## attr(,"class")
## [1] "histogram"# Using the raster object with different breaks
plot(tempRast,
breaks = c(-40,-20,0,20,40),
col = brewer.pal(4, "Spectral"))
# Reverting colours
plot(tempRast,
breaks = c(-40,-20,0,20,40),
col = rev(brewer.pal(4, "Spectral")))
1.4.2 Challenge!
Using the annual mean precipitation raster, create a plot that customises the color map with 5 breaks.
1.4.3 Cropping Rasters
## Croping rasters based on extent
usa_box<-c(-124.7258, -66.94989, 24.49813, 49.38436)
plot(tempRast)
#plot temperature raster
plot(tempRast)
#Define the extent of the crop by clicking on the plot
# cropbox1 <- drawExtent()
cropbox1<-c(-124.7258, -66.94989, 24.49813, 49.38436)
#crop the raster, then plot the new cropped raster
Tempcrop1 <- crop(tempRast, cropbox1)
Precicrop1 <- crop(precipRast, cropbox1)
#Plot the cropped extent
plot(Tempcrop1)
## Layering rasters
plot(Tempcrop1,
col=grey(1:100/100), #create a color ramp of grey colors
legend=F)
plot(Precicrop1,
breaks = c(0,1000,2000,3000,4000),
col = rev(brewer.pal(4, "Spectral")),
alpha=0.4,
add=T)
1.5 Challenge!
Create a precipitation and temperature plot of Australia
1.6 Raster calculations
#Calculate NPP using the Miami model.
npp.mia.t <- 3000/(1+exp(1.315-0.119*tempRast))
npp.mia.p <- 3000*(1-exp(-0.000664*precipRast))
nppRast<-min(stack(npp.mia.t,npp.mia.p))
plot(nppRast, main="NPP")